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ABSTRACT 

Most of super-Earths detected by the radial velocity (RV) method have sig- 
nihcantly smaller eccentricities than the eccentricities corresponding to velocity 
dispersion equal to their surface escape velocity (“escape eccentricities”). If or¬ 
bital instability followed by giant impacts among protoplanets that have migrated 
from outer region is considered, it is usually considered that eccentricities of the 
merged bodies become comparable to those of orbital crossing bodies, which are 
excited up to their escape eccentricities by close scattering. However, the eccen¬ 
tricity evolution in the in situ accretion model has not been studied in detail. 
Here, we investigate the eccentricity evolution through A-body simulations. We 
have found that the merged planets tend to have much smaller eccentricities 
than the escape eccentricities due to very efficient collision damping. If the pro¬ 
toplanet orbits are initially well separated and their eccentricities are securely 
increased, an inner protoplanet collides at its apocenter with an outer proto¬ 
planet at its pericenter. The eccentricity of the merged body is the smallest for 
such conhguration. Orbital inclinations are also damped by this mechanism and 
planets tend to share a same orbital plane, which is consistent with Kepler data. 
Such efficient collision damping is not found when we start calculations from 
densely packed orbits of the protoplanets. If the protoplanets are initially in the 
mean-motion resonances, which corresponds to well separated orbits, the in situ 
accretion model well reproduces the features of eccentricities and inclinations of 
multiple super-Earths/Earth systems discovered by RV and Kepler surveys. 

Subject headings: planets and satellites: dynamical evolution and stability - planets 
and satellites: formation 
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INTRODUCTION 


About 60 close-in super-Earths (< SOMq, Mt^ is the mass of the Earth) have been 
discovered by the radial velocity method so farl^. Fabrycky et ah (2014) showed that the 
Kepler survey found more than 818 super-Earth-sized (< OR®, R® is the radius of the 
Earth) candidates in 333 multiple systems. These planets can be formed either by 1) type 
I migration of the full-sized planets that have formed in outer regions (e.g., Cresswell & 
Nelson 2006, Cossou et al., 2014), 2) in situ accretion of planetesimals that formed there 
or have migrated from outer regions due to aerodynamical gas drag (e.g., Raymond et ah, 
2007; Chiang & Laughlin 2013), 3) in situ accretion of protoplanets that have migrated from 
outer regions due to type I migration (Terquem & Papaloizou 2007; Ogihara & Ida 2009; 
Ida & Lin 2010). Model 1 has a difficulty of why the full-sized planets were able to avoid 
runaway gas accretion. In models 2 and 3, if the growth beyond a critical core mass occurs 
after disk gas depletion, the runaway gas accretion is avoided, while observed pile-ups near 
mean-motion resonances are not easy to be explained. In model 2, accumulation of large 
amount of planetesimals may be a difficultjQ. In model 3, type I migration may be able 
to bring larger amount of solid materials to inner regions, as explained below, although 
the total mass of predicted super-Earths may still be smaller than the observed one (e.g., 
Ida & Lin 2010, Ida et al., 2013). In model 3, giant impacts among protoplanets that have 
migrated from outer regions occur after disk depletion (giant impacts could also occur in 
model 2). It is often considered that orbital eccentricities resulted in by giant impacts are 


^http://exoplanets.org 

^ Note that Chiang & Laughlin (2013) assumed that the feeding zone width is as large 
as orbital radius itself, while it is usually set to be several to ten Hill radii. So, much larger 
planetesimal surface density would actually be required than that estimated by Chiang & 
Laughlin (2013). 
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larger than the observed values (see below), which is also a problem for model 3. In the 
present paper, however, we will show that eccentricities resulted in by giant impacts should 
be as small as the observed level. 

The details of model 3 are as follows. Type I migration is halted at the disk inner edge 
and subsequently migrating protoplanets could be trapped by mean-motion resonances of 
preceding one stopped at the disk edge. Terquem & Papalouzou (2007) showed that type 
I migration is too fast for protoplanets to be captured by resonances at hrst encounters. 
They are trapped in resonances after close scattering and coagulation near the disk edge, 
resulting in a few coagulated planets in resonant orbits. Their orbits are stable even after 
disk depletion. It is inconsistent with data of Kepler candidates that most of multiple 
systems are off-resonant, unless additional effects to deviate the planets from the resonant 
conhgurations are applied (e.g., Papaloizou 2011). Model 1 also requires a similar process 
to realize off-resonant orbits (Goldreich & Schlichting 2014). 

Ogihara & Ida (2009) found that if type I migration rate is considerably reduced 
from that derived by Tanaka et al. (2002), protoplanets are resonantly trapped at first 
encounters, resulting in resonant systems consisting of a large number of protoplanets. In 
the presence of disk gas in which the eccentricities are strongly damped by planet-disk 
interactions (e.g., Tanaka & Ward 2004), the protoplanets’ orbits are stable. However, if the 
number of resonant protoplanets exceeds a critical value, the systems can become unstable 
after disk gas depletion (Matsumoto et al. 2012). The following orbit crossing and giant 
impacts result in off-resonant multiple super-Earth systems (Ogihara & Ida 2009), which 
could be consistent with the off-resonant Kepler systems. 

The observed orbital eccentricities and inclinations constrain the formation model 
as well as semimajor axis distributions. Eccentricities are estimated by radial velocity 
measurements, while mutual inclinations are constrained in multiple systems by transit 
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detection. In model 1, the planetary orbits should be almost coplanar and circular. Higher 
eccentricities and inclinations are expected in models 2 and 3. In this paper, we will discuss 
the latter model in comparison with the observation, because model 1 has a difficulty of 
runaway gas accretion and we will show that the eccentricities and inclinations in the latter 
models are not actually high. 


While some gaseous giant planets have eccentricities as large as 0.9, super-Earths and 
Neptune-type planets generally have smaller eccentricities than gas giants (Mayor et al. 
2011). Scattering between similar-sized planets can excite their velocity dispersion up to 
their surface escape velocities (uesc) during their assemblage stages in a gas free conditions 
(e.g., Safronov 1969; Aarseth et ah, 1993; Kokubo & Ida 2002). If the velocity dispersion 
exceeds Uesc, the collision cross section becomes larger than the scattering cross section, so 
that the excitation of eccentricities is saturated at Vesc- The corresponding eccentricity is 
given by Cesc ~ Vesc/vK where vk is the Kepler velocity, which we call “escape eccentricity”: 




2(Affc + Ml) a 


M* 


Rk Ri 
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where M* is the mass of the central star, Mq is the solar mass, Mk and Mi are masses of 
protoplanets, M^ is the Earth mass, Rk and Ri are physical radii of protoplanets, p is the 
material density of protoplanets, and a is the semimajor axis of protoplanets (Kokubo & 
Ida 2002). 


However, most of observed super-Earths have smaller eccentricities than their escape 
eccentricities. In Figure [H we show observed eccentricities of planets observed by RV 
method scaled by their escape eccentricities, where we omit planets inside of 0.1 AU because 
their eccentricities can be damped by tidal interactions with their host stars (Goldreich & 
Soter 1966). Since planetary radius is needed to estimate Cesc (RV observations give only 








6 


planetary minimum mass), we assume a mass-radius relationship for planets whose densities 
are not known. Figure [T] shows eje^sc of observed 22 planets as a function of planetary 
mass. Since the mass-radius relation has large uncertainty, we tested three models. In panel 

(a) , densities are given by 3 g/cm^, independent on planetary mass. This is equivalent 
to the mass(Mp)-radius(i?p) relation, Mp/Mf^ = O.54(i?p/i?0)^, where i?® is the Earth 
radius. This simple model is often used in N-body simulations for rocky planetesimals or 
protoplanets (e.g., Kokubo et ah, 2006). This panel shows that super-Earths have smaller 
eccentricities than their escape eccentricities. However, it is observationally suggested that 
larger super-Earths tend to have lower bulk densities. Lissauer et ah (2011) derived an 
empirical mass-radius relation as Mp/M® = (i?p/i?®)^'°® by htting the planet in the solar 
system. With this mass-radius relation, the density of SOM® planet is 1.2 g/cm^. In panel 

(b) , this relation is used. Even with different mass-radius relations, the e/cesc distribution 
in panel (b) is similar to that in panel (a). Since the escape eccentricity is proportional to 

the different density does not signihcantly affect the result. From the lower density 
of observed larger exoplanets, Wu & Lithwick (2013) and Weiss & Marcy (2014) derived 
mass-radius relations with stronger density-dependence on mass. Wu & Lithwick (2013) 
derived Mp/M® = 3(i?p/i?®) for 1.6 < Rp/R^ < 7 planets (Weiss & Marcy (2014) derived a 
similar relation). This model is adopted in panel (c). Although 4 planets have e > Cesc, they 
can be e < Cesc within the error bars. Furthermore, all of them have Mp > 18M® and their 
densities are estimated to be p < 0.4 g/cm^, which may be lower-estimated. Thereby, We 
conclude that the eccentricities of super-Earths are e < Cesc • This means that super-Earths 
were formed through in situ coalescence followed by some eccentricity damping or not 
formed through the in situ coalescence. 

Planetesimal accretion near 1 AU has been extensively studied. Planetesimals grow 
through runaway growth in early stage (e.g., Wetherill & Stewart 1989; Kokubo & Ida 1996) 
and the oligarchic growth follows (e.g., Kokubo & Ida 1998). In these stages, protoplanets 
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grow up, accreting surrounding planetesimals crossing their orbits from various directions. 
Planetesimals prevent protoplanets from orbital crossings between protoplanets thanks to 
the dynamical friction. This process is referred to as the orbital repulsion (Kokubo & Ida 
1995). After most protoplanets accrete planetesimals in their feeding zones, the dynamical 
friction becomes no more available. After disk gas is depleted, the planet-disk interaction 
is not available for the eccentricity damping as well. Then, orbital instability occurs. This 
stage is called as the giant impact stage. The accretions of protoplanets in the giant impact 
stage are investigated in several papers (e.g.. Chambers & Wetherill 1998; Agnor et ah, 
1999; Raymond et ah, 2004; Kokubo et ah 2006). They are successful in reproducing 
terrestrial planets in the solar system in some aspects, such that the Earth and Venus mass 
planets are formed around 1 AU, if the range of initial semimajor axes of protoplanets is 
relatively restricted (Hansen 2009). 

The eccentricities of the formed Earth-mass planets are usually e ~ 0.1 (Kokubo et 
ah, 2006). While they are larger than the current free eccentricities of Venus and the 
Earth ~ 0.02 — 0.03, they are 3 times smaller than Cesc- Although some external damping 
mechanisms, such as dynamical friction from residual planetesimals (O’Brien et ah 2006) 
or planet-disk interaction (dynamical friction from disk gas; Kominami & Ida 2004) need to 
be taken into account to reproduce the current eccentricities of Venus and the Earth, the 
low eccentricities of observed super-Earths can be potentially explained by the accretion 
through the giant impacts in gas-free environment. We will make clear why eccentricities of 
formed planets become smaller than Cesc- 

For transit of multiple planets at a < Oout to be detected, mutual inclinations must be 
within 2.6(R*/i?0)(a/O.l AU)“^ degree. Fabrycky et al. (2014) suggested that the typical 
mutual inclination of multiple super-Earths in Kepler candidates lies hrmly in the range 
i = 1.0 — 2.2 degree. On the other hand, numerical simulations show that close scattering 
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between planetesimals realize e ~ (e.g., Ida & Makino 1992; Shiidsuka & Ida 1999). For 

e ~ Cesc, i ~ 4sc ~ eesc/2, which is 5.4 degree for a lOM® super-Earth at O.lAU. This is 
larger than the mutual inclinations of Kepler systems. However, N-body simulations show 
that Earth mass planets formed around 1 AU normally have i ~ 3 degree (e.g., Kokubo et 
al. 2006), while 4sc ~ 8.6 degree. This result suggests that small inclinations of the Kepler 
systems can be realized by the in situ accretion in gas-free environment. 

In this paper, we study how the hnal velocity dispersion becomes smaller than the 
escape velocity to account for the small eccentricities of observed super-Earths through 
A^-body simulations of in situ accretion of planets in giant impact phase near the central 
star and analytical arguments. In Section |2l we outline the numerical methods and initial 
conditions of protoplanets. Our results of A^-body simulations are presented in Section 
131 From 35 simulations in section 13.11 we hnd the eccentricity and inclination of the 
largest planets are typically much smaller than the escape velocity. Because planetesimals 
and planet-disk interaction are not included in the simulation, the low eccentricities and 
inclinations are not due to dynamical friction. In section 13.21 we find collisions in the giant 
impact stage tend to occur when the angles between pericenters are around 180 degree and 
they are responsible for the low eccentricities and inclinations. We summarize the results in 
Section 01 


2. NUMERICAL MODEL 


We perform A-body simulations of the planet accretion starting from protoplanets 
without small planetesimals. Orbital evolution of protoplanets is obtained by the numerical 
integration of 


d'^rj 

dt"^ 
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where G is the gravitational constant, r* are the coordinate of the i-th. protoplanets 


numbered from the innermost, and Vij is the relative distance of the planets i and j. In our 


calculations, the central star has a solar mass. The numerical scheme is the fourth-order 
Hermite scheme. We assume perfect accretion, i.e., planets always accrete without bouncing 
nor fragmenting, and the total momentum is conserved at a collision. Every time a collision 
is taken place, we record the orbital elements and masses of protoplanets before and after 
the collision. 

We perform 2 sets of simulations: 5 — 16 non-equal-mass protoplanets (which is 
hereafter referred to as “iV-body set”) and three equal-mass protoplanets (“three-planet 
set”). In a standard case (case A; 20 runs) of the A^-body set, we distribute protoplanets in 
a range from 0.05 AU (= Oi) to 0.29 AU. The number of the planets N = 16. Their total 
mass is Mtot = 17.3M0. Individual masses are given by 



(3) 


with El = 100 gcm“^. Although these protoplanets may have migrated from outer regions, 
we used a formula for isolation masses in in situ oligarchic growth (Kokubo & Ida 2002). 
Orbital separations (6) are lOrn where rn is the Hill radius. The initial individual masses 
do not affect the results as long as 6 ~ lOrn. We also performed similar simulations 
with different Si and N (accordingly, Mtot is also different) with the same b: case B, 

C, and D (5 runs for each set). Planetary physical radii are calculated using a material 
density of p = 3 gcm“^. The initial eccentricities and inclinations of protoplanets are 
given by the Rayleigh distribution. The dispersions of eccentricity and inclination are 


(g 2 )i /2 = 3 X 10 ^(Si/100 gem ^)^/^ and = 6 x 10 “^(Si/lOO gem ^)^/^ radian. 


The initial conditions are summarized in Table [T] The simulations follow the evolution of 
protoplanet systems for 10® Kepler time of the innermost planet. In some simulations, we 


calculate 3 x 10® Kepler time and confirm that the resultant planet are stable. 
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In the next set of simulations, we perform many runs using equal-mass three-planets 
for statistical surveys. These simple-settings enable us to control planetary masses and 
semimajor axes of colliding bodies. In the three-planet cases, we give the semimajor axes 
of middle planets ( 02 ). Inner planets and outer planets are set at 01,3 = 02 ± hr-n. The 
Hill radius is given by rn = ( 2 Mp/ 3 M*)^/^a 2 ~ 1.26 x lO“^(Mp/M 0 )^/^(a 2 /O.l AU) AU. 
The planetary radius (i?p), planetary mass (Mp), the semimajor axis of the middle planet 
( 02 ), the orbital separations normalized by the Hill radius ( 6 ), and the initial eccentricities 
(cini) are free parameters. The orbits of planets are coplanar. We perform 11 cases in total, 
and we calculate 100 runs in each case changing initial orbital angles of the protoplanets 
randomly. The initial conditions are summarized in Table [2l We also calculate systems 
composed by non-zero inclination planets (a 2 hni < 2 rH), and confirm that collisions between 
inclined planets show the same tendency of collisions between planets in coplanar orbits, 
although the results in non-coplanar cases are not presented in this paper. 


3. RESULTS 

We hrst present the results of the A^-body set {N = b — 16). We focus on the 
eccentricities and longitudes of pericenter before and after collisions to investigate 
eccentricity evolution through collisions. Next, we show the results of the three-planet set, 
to study the dependences of eccentricity evolution on initial conditions. Through these 
calculations, we explain intrinsic dynamics to cause the efficient collisional damping for the 
systems starting from moderately separated orbits. 
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3.1. Results of iV-body Set 

The typical orbital evolution is shown in the left penal of Figure |2l This hgure shows 
time evolution of semimajor axes, pericenters, and apocenters of planets for 1.0 x 10^ yr. 
In this calculation, six planets are hnally formed. They have hnal masses between l.SM® 
and 4.5M0. Their escape eccentricities are Cesc = 0.079 - 0.23 (eqnation ([1])). The hnal 
eccentricities of planets are between 0.015 and 0.047, which are mnch smaller than Cesc- 
The (e)/eesc of the largest, the second largest, and other planets in case A, B, C and D are 
summarized in Table [3l In all cases, (e)/eesc is less than nnity. In particular, for the largest 
bodies, (e)/eesc is only 0.1-0.2 except case D in which Mtot is extremely small (~ O.I 5 M 0 ). 

These e/cesc < 1 featnres are caused only by collisions. The middle panel of Figure 
[2] is the closeup of the eccentricity evolution of protoplanets at a collision. The fourth 
innermost planet with 2 .IM 0 and the hfth one with I. 2 M 0 collide at t ~ 4.6 x 10^ yr. The 
eccentricities of the inner and onter planets are 0.066 and 0.10 just before the collision, 
which are comparable to Cesc ~ 0.1. However, the eccentricity of the merged body is 0.012 
just after the collision, which is an order of magnitude smaller than Cesc- 

We show the orbits of two planets just before the collision in the right panel of Figure 
[21 The locations of the collision and their pericenters at the collision are shown by hlled 
circles and crosses, respectively. The azimnthal velocity at the apocenter of the inner planet 
is given by a/G'M*( 1 — ei)/ai(l -1- ei) ~ uk(1 — ei), where uk is the Keplerian velocity at 
the collision location, while that at the pericenter of the outer planet is uk(1 + 62 ). Since 
they have similar masses and eccentricities, the velocity of the merged body should be 
~ uk due to conservation of total momentnm, which means that the orbit of the merged 
body is nearly circnlar. If the orbital separation is comparable to radial excursion dne to 
the eccentricities, collisions occur only when the apocenter of the inner planet meets the 
pericenter of the outer planet (if their pericenters are aligned, their orbits never cross). 
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For such orbital separation, the collisional damping for eccentricity is always very efficient. 
If the orbits of protoplanets are nearly circular and well separated, the eccentricities are 
excited only by secular perturbations. The eccentricities are secularly increased until the 
apocenter distance of the inner planet approaches the pericenter distance of the outer planet 
and a collision between them occurs. 

More detailed analysis can be done using Laplace-Runge-Lenz vector. The mass- 
weighted total Lenz vector is conserved during scattering and even collisions under Hill’s 
approximation (Nakazawa & Ida 1988). According to the conservation, when two planets k 
and I collide and are merged into a planet kl, the eccentricity of the merged body (e^;) is 
written as 


(Mfc + Mifeh = Mlel -f- ef 2MkMiekei cos {wk - ru/), (4) 

where w are longitudes of pericenters of the bodies. When longitudes of pericenters of 
bodies are randomly distributed, the average of eh becomes 

{Mk + M,feh = Mlel + Mfel (5) 

This equation means that the eccentricity of the merged body is comparable with those of 
the colliding bodies; in the case of Mk = Mi and Ck = ei, Cki = ekj^2. The approximation 
of random w is valid if we consider a radially packed distribution of bodies. The validity 
of equation ([5]) is confirmed by the A^-body simulations for random velocity evolution of 
packed planetesimals (Ohtsuki 1992). 

However, as already mentioned, in the case of giant impacts of protoplanets that have 
initially well separated orbits, Acu = Wk — voi may be ~ 180 degree. In this case, equation 
(jH) implies 


(A4 + ~ M'lel + M'lef - 2MtMiekei. 


( 6 ) 
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When Mk ~ Mi, eh ~ {ck — e«)^/4, which is much smaller than eh ~ (e^ + ef)/4 given by 
equation (E]). 

For the collision in Figure [21 the inner planet has Mk = 2 .IM 0 and = 0.066 and 
the outer planet has Mi = I. 2 M 0 and ei = 0.10. We found Atxj = 171 degree. If we use 
equation (E]), the estimated value of is ~ 0.056. However, it is ~ 0.0080 with equation 
(14|) . which is much more consistent with the orbital integration. 

The collisions like the right panel of Figure |2] occur commonly between separated 
protoplanets. In case A, there are 203 collisions in 20 runs. The panel A of Figure ID shows 
the distribution of Aw obtained in the 203 collisions. It clearly shows that collisions tend to 
occur around Aw = 180 degree. Despite the difference in the number, masses, semimajor 
axes of protoplanets, similar peaks at Aw = 180 degree are found in the other cases. The 
numbers of collisions are 13 in case B, 24 in case C, and 21 in case D. The mean values and 
variances of Aw are 180 ± 53 degree in case A, 177 ± 30 degree in case B, 178 ± 63 degree 
in case C, and 187 ± 55 degree in case D. 

This means that the eccentricity of merged protoplanets tend to be much smaller 
than the estimation in equation ([5]). This efficient eccentricity damping was mentioned in 
Raymond et ah (2006), although they did not analyzed the concentration of Aw around 
180 degree. 

If orbit crossing still continues, the damped eccentricity is excited again to ~ Cesc- 
However, because timescales (tctoss) for orbital instability to start sensitively depend on the 
initial orbital separations (e.g.. Chambers et ah 1996), Tcross of a system can jump up by 
several orders of magnitude at a collisional merging (see Figure 3 in Ida & Lin 2010). After 
that, the system becomes stable during main-sequence lifetime of the host stars and the 
damped eccentricities are remained. 
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We find inclinations are also significantly damped throngh collisions. Figure [5] shows 
the inclinations just before and after collisions. The velocity component normal to the 
invariant plane depends on the ascending node. When a collision occurs at the ascending 
node and descending node of colliding bodies, the velocity component normal to the 
invariant plane of the merged body is much smaller than those of colliding bodies. 

The largest planets formed in case A have the averaged inclination of (i) = 1.2 ± 1.8 
degree. Since the largest planets have a mean mass ~ 4.4M0 and semimajor axis ~ 0.2 AU 
(Table [3]), iesc ~ eesc/2 — 6 degree. The inclinations of the largest planets are considerably 
smaller than igsc- In other A^-body set cases, formed planets also have smaller inclinations 
than iesc- The inclinations of the largest bodies are {i) = 4.7 ± 1.9 degree, 1.9 ± 3.0 degree, 
and 0.70 ± 0.23 degree in case B, case C, and case D, respectively. Mutual inclinations in 
case A and C agree with those of Kepler planets, i = 1.0 — 2.2 degree (Fabrycky et ah 
2014). Because of the inclination damping, the hnal planetary systems tend to be coplanar. 
The means and variances of the typical mutual inclinations among all planets in a system 
are (irei) = 1-3 ± 1.7 degree, 4.4 ± 7.0 degree, 0.78 ± 1.3 degree, and 6.8 x 10“^ ± 8.6 x 10“^ 
degree in case A, case B, case C, and case D. Because of larger semimajor axes of planets 
in case B, iesc is larger. Accordingly, (i) of the largest planets and (irei) are larger in case 
B than those in the other cases, although (i) and (irei) are still < iesc- In general, angular 
momentum dehcits (AMDs) are increased from initial values by scattering. However, the 
increase is not so signihcant except in case B (Figure [3]). 

In case A, C, and D, protoplanets tend to collide with the neighboring protoplanets 
rather than undergo global orbital instability. In the proximity of their host stars, the ratio 
of Hill radii to physical radii is small, so that scatterings are less dominated over collisions 
than in outer regions. In other words, Cesc and iesc are smaller for smaller semimajor axis. 
Furthermore, e and i are signihcantly smaller than Cesc and iesc- Therefore, e and i can be 
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very small through collisional damping in close-in regions. Our results are not affected by 
the assumption of perfect accretion. It was shown that the eccentricities in hybrid N-body 
and SPH simulations allowing collisional fragmentation by Kokubo Sz Genda (2010) do not 
differ from those obtained in perfect accretion simulations. 


3.2. Results of Three-Planet Set 

The very effective eccentricity damping comes from the concentration of Acn on 180 
degree in collisions. Let = lAtn — 7r|. Assuming 1, equation (jl]) becomes 

(Mfc -I- Mife\i = (MfcCfc — Mieif -|- MkMiCkeie^, (7) 

when Mk = Mi and Cfc = e; ~ Cesc, 

^kl /0\ 

- ~ — ~ ---, (^oj 

Cesc 2 115 degree 

If Cro ^ 10 degree, the eccentricity of the merged body is an order of magnitude smaller 
than Cesc- 

In this section, we investigate how the concentration occurs through simpler three- 
planet simulations. We change planetary radii, planetary masses, initial semimajor axes 
of the middle planets, initial orbital separations, and initial eccentricities. The initial 
conditions are summarized in Table |2j The general features of results of three-planet 
calculations are basically similar to those of A-body set in the previous section. The 
Azu distributions are peaked at 180 degree, eki is significantly smaller than Cesc, and the 
estimated eccentricities in equation (|1]) agree with e^i- 

With initial spacing we use, the system in the three-planet set readily becomes unstable 
and orbit crossing starts. In this case, the system enters stable state after a hrst collision 
and the eccentricity damped at the collision is not usually excited any more. Thereby, Azu, 
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the degree of eccentricity damping and their dependences on initial conditions are better 
described than in iV-body set. 

Case 3A is the standard case of three-planet calcnlations. In this case, the second 
innermost planet is at 0.1 AU. All planets have the same masses (IM®), radii (li?®), and 
initially circular orbits (e = 0). The initial orbital separations are given as 4rH, which is 
equal to 5.04 x 10“^ AU. The distribution of differences between pericenters of colliding 
bodies {Avu) in 100 simulations is shown in Figure |6l The Azu distribution in case 3A is 
peaked at 180 degree, in the same manner as case A (Figure llj). The variance of the Aw 
distribution (a^) is 17 degree. 

Typical evolution of eccentricities and arguments of pericenters in case 3A is as follows. 
Eccentricities of planets increase by mutual scatterings, and collisions occur not long after 
their eccentricities exceed the eccentricities required for orbital crossing (ccross = da/2a). 
When the orbits of two planets hrst become able to collide with each other, their orbits 
should have e ~ Ccross- However, the planets usually undergo close encounters and their 
e are excited from Ccross before an actual collision. In case 3A, the mean of the larger 
eccentricity between colliding planets before the collision is 0.094, which is larger than 
ecross(— 0.028), and comparable to eesc(— 0.11). The detailed eccentricity evolution shows 
that colliding planets, which have typically e = 0.041 ~ l.Sccross are pumped up to above 
eccentricity just before the collision. The mean eccentricity after the collision is reduced to 
0.015. Substituting Aw = 180 ± 17 degree, = Mi, and ~ ei into equation (jl]), we get 
^ki/^k = 0.15, which agrees well with the numerical value 0.015/0.094 = 0.16. 

In the following, we discuss Aw distribution using the pericenter dispersion First, 
we analytically estimate If eccentricities secularly increase from zero, the collision 
between two planets becomes possible when the apocenter of the inner planet (Qi) contacts 
with a pericenter of the outer planet (^ 2 ) with Aw = 180 degree. In the following analysis. 
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we neglect eccentricity excitation from Cesc for simplicity. Althongh the assnmption is not 
relevant enongh, the analytical discnssion neglecting the excitation well reprodnces the 
numerical results. 


If we take into account physical radii of planets, the collisional point can rotate by an 
angle from the pericenter of the outer planet (Figure [7]). In this case, 


02(1 


R 2 — ®i(l + Cl) + i?i. 


(9) 


1 + 62 COS 

where Ri and R 2 are planetary physical radii of the inner and outer planets, respectively, 
and the true anomaly of the outer planet is equal to e^. Under the assumption that planets 
collide at Qi, the angle + Azu = vr. When we assume ~ g 2 3> Rtot = Ri + R 2 , 

1 + 62 Rtot 


62 COS = 62- 


1 — 62 02 


Assuming 62 = 6 1, 1 radian) is 


f 


2Rtot 

602 


( 10 ) 


( 11 ) 


With 6 ~ 6 cross = do/ 2 o 2 = 6 rH/ 2 o 2 , equation flTTll reads as 


2i ' 


br 


H 


= 15|A 


1/2 


Mrr, 


Mo 


1/6 


0.1 AU 


-1/2 


( 12 ) 


[degree]. (13) 


The estimation of agrees well with cr^ in case 3A (<7^=17 degree). For collisions with 
e > 6cross, simulated 0 ^ is larger than given by equation (IT^ . resulting in inefficient 
eccentricity damping. In the case of a dense orbital distribution of planetesimals. Aw is 
uniformly distributed without any concentration at vr. Then, the collisional damping is not 
effective and the eccentricities of merged bodies are similar to those during orbit crossing. 

To check the validity of equation ([13]), we perform additional runs. Figure [8] shows the 
results of different i?p. We found =17, 24, and 39 degree in case 3A (i?p = i?®), case 3B 
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{Rp = 10^/^i?®), and case 3C (-Rp = lOi?®), respectively. They agree with corresponding 
estimations, ~ 15, 27, and 47 degree. We also perform calculations with changing 
planetary radii and masses in case 3D and case 3E, keeping R^/da constant, which means 
that is constant (Figure |9]). The resultant is 16 — 18 degree in all cases, while Rp 
changes a factor of 10. These results also indicate that Mp affects through da = brn, 
and (Tro is proportional to Mp^^^. 

Next, we change da by changing 02 and b with hxed Mp. The dependence on 02 is 
shown in Figure HD) For 02 = 0.1^/^ AU (case 3F), 1 AU (case 3A) and O.l^A AU (case 3G), 
cr^ obtained by simulations are 22, 18 and 17 degree. The estimated by equation flT^ is 
not relevant enough in case 3G, because eccentricities are more highly pumped up before 
collisions. Figure ITT] shows the results with 6 = 4 (case 3A), 5 (case 3H) and 6 (case 31). 
Although the timescale for orbital instability to start is very different (e.g., Ghambers et ah 
1996), cr^ are similar: 17 degree in case 3A (6 = 4), 16 degree in case 3H (6 = 5), 21 degree 
in case 31 (6 = 6), which are consistent with ~ 15 ( 6 / 4 )“^A degree (equation ffT3|) h These 
results show that concentrations of Aw at vr with (7^ — 20 degree is quite common as long 
as close-in regions (< 0.3 AU) are considered. 

The above results are applied to the systems in which mean orbital separations are 
larger than 2\/3rH and orbital crossing doe not occur until eccentricities are gradually 
increased by distant perturbation. In a system of equal-mass bodies with surface density S 
at a, the mean orbital separation is 



(14) 


If we consider early stages in which the systems consist of plenty of small planetesimals, 6 


is far smaller than 2\/3. Then, the concentration of Azu does not occur and the collisional 


damping should be weak. If e is set such that the radial excursions are larger than orbital 
separations (ea 2 > 6rH), a situation is similar. We set 002 > 6rH in case 3K. As expected. 





we obtain a relatively large value of <7^ (=56 degree) in this case, because collisions occur 
regardless of the directions of pericenters and Azu is no longer concentrated (Figure [T2]). 


4. SUMMARY AND DISCUSSION 

We have investigated the eccentricity damping through the giant impacts of the 
protoplanets in the proximity of the host stars. First, we performed 20 runs of A-body 
simulations of protoplanets starting from 16 bodies of about Earth-mass at 0.05 — 0.29 AU 
with orbital separation of 10 Hill radii (“A-body set”). We also performed simulations 
of runs with more and less massive bodies and runs at larger semimajor axes. We have 
conhrmed that eccentricities of formed planets are signihcantly lower than the eccentricity 
corresponding to velocity dispersion of their surface escape velocity (“escape” eccentricities 
Cesc; see Section IXT]) . For an Earth-mass body at O.lAU, Cesc ~ 0.1. During orbital crossing, 
eccentricities increase due to the mutual scatterings among protoplanets and reaches ~ Cesc- 
However, the eccentricities are damped by an order of magnitude at a collision. When 
the orbits of the protoplanets are relatively separated and their eccentricities are secularly 
increased, the differences between the pericenters of colliding planets tend to be Azu ~ 180 
degree, i.e., the collisions occur at the apocenter of the inner body and the pericenter of 
the outer body. Since the azimuthal velocity of the inner body is slower than the local 
Keplerian velocity and that of the outer body is faster and the two bodies have similar 
masses, the velocity of the merged body should be close to the local Keplerian velocity, 
which means that the orbit of the merged body is nearly circular. We also described more 
detailed discussion on why the collision damping is so efficient, using conservation of Lenz 
vector in Hill’s approximation. 

The damped eccentricities are again excited up to ~ Cesc if orbital crossing continues. 
However, after some merging, the planets become isolated from one another and the planets 
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are remained in stable orbits with e <C Cesc- We found that the inclinations of protoplanets 
are also damped through collisions. The mutual inclinations among formed planets in the 
massive systems in close-in regions (case A and C) are (i) = 1.3 ± 1.7 degree in case A 
and {i) = 0.78 ± 1.3 degree in case C, respectively, which agree with those of observed 
super-Earths, i = 1.0 — 2.2 degree (Fabrycky et al. 2014). 

Next, we performed three-planet calculations (“Three-planet set”). With these simple 
systems, more detailed analysis on the collision damping can be done. The eccentricity 
of a merged body is given approximately by e ~ (ct^/llh degree)eesc where is the 
width of the concentration of Azu around 180 degree (equation ([8])). If the concentration 
is high (e^ is small), e can be much smaller than Cesc- Through the analytical argument 
and the results of the three-planet calculations, we found that is approximated as 
Cro ~ 2(2i?p/(ia)^/^ where Rp is planetary physical radius and da is initial orbital separation. 
That is, e ~ 0.13{Rp/R(^Y^‘^{da/0.05 AU)“^/^eesc- Note that in runaway and oligarchic 
stage, the collision damping is weak and e ~ Cesc, because da is very small and accordingly 
Azu is uniformly distributed. 

Volk & Gladman (2015) studied orbital stability of multiple super-Earth systems 
analogous to the systems discovered by Kepler and found that in some systems, global 
instability occurred and the collision velocities are larger than 2uesc (the collisions are 
disruptive), while in our simulations, no global instability occurred and collision velocities 
are smaller than 1.5uesc- Although Volk & Gladman (2015) did not explain the conditions 
for global instability, initial planetary mass distribution may have caused the difference. In 
their calculations, the maximum mass ratios among initial planets are generally much larger 
than ours. If a small-mass body collides with another small-mass one after scattering by a 
massive planet, the collision velocity can be larger than Uesc of the small bodies. While the 
results of Volk & Gladman (2015) can be applied to some of Kepler systems, our results can 
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be applied to other systems. Detailed investigation of scattering in systems with initially 
large mass ratios is left for future work. 

Figure [T3j4 shows the eccentricities e and masses Mp of the formed planets in iV-body 
set (cases A and C). Observed data of super-Earths in Figure [T] (a) are also plotted for 
comparison fFigure [T^h We only plot planets at a > 0.1 AU, because eccentricities of 
planets at a < 0.1 AU may be damped by tidal dissipation. The eccentricities of the formed 
super-Earths in our calculations are e < O.Scesc except a few planets. In the observed 
data, when we adopt p = 3 g/cm^, all planets have e < Cesc- The data are consistent with 
our result, although our results show slightly smaller e. Even if we use a more realistic 
mass-radius relation based on the Solar system planets derived by Lissauer et al. (2011), 
the result is hardly changed. With the mass-radius relation of Wu & Lithwick (2013) that 
produces lower bulk density, e > Cesc for planets with > 2OM0. However, p ~ 0.4 g/cw? for 
these planets may be too small. Furthermore, the best-£t eccentricities of planets detected 
with the low signal-to-noise ratio and the small number of observations by RV surveys 
tend to be larger than the true values (Shen & Turner 2008). Therefore, we conclude that 
the eccentricities of observed close-in super-Earths are not inconsistent with our results. 
Note that our results with p = 3 g/cm^ are equivalent to the results with other density 
and semimajor axis according to a scaling. The dynamical process is scaled by the ratio 
of the geometrical cross section to the scattering cross section with radius th, which is 
proportional to p~‘^/^a~‘^. So, our results at a = 0.5 AU correspond to the case of p = 1 
g/cm^ and a = 0.035 AU. For example, GJ 667C c has e = O.OTcesc and 4.2M®. GJ 667G is 
a member of a triple stellar system. The eccentricity of GJ 677G c could be affected by GJ 
667A and GJ 667B. Therefore, the eccentricities of super-Earths in observed data may be 
consistent with those obtained in our N-body simulations. 


As we described in section 1, three models have been proposed for formation of close-in 
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super-Earths: 1) type I migration of the full-sized planets that have formed in outer 
regions, 2) in situ accretion of planetesimals that formed there or have migrated from outer 
regions due to aerodynamical gas drag, and 3) in situ accretion of protoplanets that have 
migrated from outer regions due to type I migration. As explained in Introduction, none of 
these can completely explain the presently known super-Earth systems. What we argued 
in the present paper is that relatively low eccentricity found in Kepler systems is not a 
negative factor for a formation model of close-in super-Earth systems via giant impacts 
of protoplanets. With this result, model 3 might look promising. However, more detailed 
discussions on different aspects are needed to clarify the origin of super-Earth systems 
discovered by Kepler. 
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Fig. 1.— The exoplanets detected by RV method are plotted. Data was extracted from 
http://exoplanets.org. The horizontal axis is the planetary mass normalized by the Earth 
mass. The vertical axis is the eccentricities normalized by their escape eccentricities, defined 
by equation m- We adopt different density models in each panel. In panel (a), we adopt 
p = 3 g/cm^ for planets whose densities are not known. In panel (b), planetary radii are given by 
Mp/M 0 = (Rp/Rq)^'^®, following Lissauer et al. (2011). In panel (c), planetary radii are given 
by Mp/M® = 3{Rp/R^), following Wu & Lithwick (2013). The planets whose semimajor axes 
are greater than 0.1 AU are plotted. We remove planets whose eccentricities are unknown and 
assumed to be 0. The horizontal dashed line is e/cesc = 1- 
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Fig. 2.— The orbital evolution of a run of case A. Left: The time evolution of semimajor axes 
(solid lines), and pericenters and apocenters (dashed lines) of all planets. Middle : Close-up of the 
left panel at a collision. The time evolution of the eccentricities of two colliding protoplanets is 
plotted. They collide at t = 4.60 x 10^ yr. Red solid, green dashed, and blue solid lines are inner, 
outer, and merged protoplanets, respectively. Right: The face-on view of the orbits of colliding 
bodies is plotted. The orbits of the inner and outer planet are shown in the solid curve and dashed 
curve, respectively. The central star is located at the origin. Circles indicate the positions of 
protoplanets. Crosses are the locations of their pericenters. The two pericenters are located in the 
opposite direction. 
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Fig. 3.— The initial AMD and final AMD of calculations in A^-body set are plotted. Circles 
those in case A, triangles are in case B, squares are in case C, and crosses are in case D. 




fraction 


29 



Q ... I.. I III I I il li I III Mil III Mil I III III Mil ill i illil li i ill il I III! il I III I li lillliMi i il 

0 90 180 270 360 0 90 180 270 360 0 90 180 270 360 0 90 180 270 360 


ACD ACiiJ ACD ACD 


Fig. 4.— The distributions of the difference between pericenters of colliding two bodies just 
before the collisions in A^-body set are shown. Bin size is 10 degree. In case A, There are 203 
collisions in total. Total numbers of collisions are 13-24, in the other cases. The mean and variance 
of Acu are 180 ± 53 degree in case A, 177 ± 30 degree in case B, 178 ± 63 degree in case C, and 
187 ± 55 degree in case D. 
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Fig. 5.— Inclinations just before and after collisions in A^-body cases are plotted in units 
of radian. The larger inclination of the colliding two bodies is denoted as ik, and iki is the 
inclination of the merged body. Circles are those in case A, triangles are in case B, squares 
are in case C, and crosses are in case D. The dashed line is iki = ik- In every collision, iki < ik holds. 
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Fig. 6.— Same as Figure HI but for case 3A, where we calculate orbital evolution of three 



planets. The mean and variance of Am are 178 ± 17 degree. 
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Fig. 7.— The orbital configuration of colliding planets is illustrated. The large filled circle 
is the central star, and the small two filled circles are planets. The two ellipses are the orbits 
of planets, and a dashed line is connecting the pericenter and apocenter of a planet. The 
angle between the pericenters is Avj. Planets collide at the apocenter of the inner planet 
(Qi), which is rotated by the pericenter of the outer planets (q' 2 )- 
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AOJ 

Fig. 8.— Same as Figure 01 but for case 3C, case 3B, and case 3A. The panels are arranged 
in order of decreasing planetary radii from the top, Rp = lOT?© in case 3C, in case 3B, 

and li?® in case 3A. The variances of Aw are 39 degree in case 3C, 24 degree in case 3B, and 17 
degree in case 3A. Decreasing Rp, Aw is more concentrated on 180 degree. 
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Atn 

Fig. 9.— Same as Figure EJ but for case 3E, case 3A, and case 3D. In these cases, we change 
planetary radii and masses keeping constant. The panels are arranged in order of 

decreasing planetary radii from the top, Rp = 10^/^i?® in case 3E, IT?® in case 3A, and 
in case 3D. The variances of Azu are 18 degree in case 3E, 17 degree in case 3A, and 16 degree in 


case 3D. 
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Fig. 10.— Same as Figure 01 but for case 3G, case 3A, and case 3F. The panels are arranged in 
order of decreasing semimajor axis from the top, 0.1^/^ AU in case 3G, 0.1 AU in case 3A, and 
02 = 0.1^/^ AU in case 3F. The variances of Aw are 18 degree in case 3G, 17 degree in case 3A, 


and 22 degree in case 3F. 
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Fig. 11.— Same as Figure IU but for case 31, case 3H, and case 3A. The panels are arranged in 
order of decreasing da from the top, da = 6rH in case 31, Sru in case 3H, and drn in case 3A. The 


variances of Aw are 21 degree in case 31, 16 degree in case 3H, and 17 degree in case 3A. 
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Fig. 12.— Same as Figure 01 but for case 3K, case 3J, and case 3A. The panels are arranged in 
order of decreasing a 2 eini from the top, a 2 eini = 5rH in case 3K, 2.5rH in case 3J, and Orn in case 


3A. The variances of Atu are 56 degree in case 3K, 25 degree in case 3J, and 17 degree in case 3A. 
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Fig. 13.— (A): Eccentricities of planets formed in our A^-body calculations. The formed planets 
that have a > 0.1 AU in the results of case A and case C are plotted to compare to Figure [H 
The filled triangles are the planets in case A, and the open triangles are those in case C. (a): 
Eccentricities normalized by their escape eccentricities with p = 3 g/crv? of observed planets are 
plotted (Same as Figure [Da). 











Table 1. Initial conditions of iV-body simulations 



N 

Si 

[gcm“^] 

ai 

[AU] 


(e2)V2 

Rl/ i?® 

^ of run 

case A 

16 

100 

0.05 

17.3 

3.16 X 10-2 

0.99 

20 

case B 

5 

100 

0.68 

25.3 

3.16 X 10-2 

1.90 

5 

case C 

8 

300 

0.05 

37.2 

5.48 X 10-2 

1.72 

5 

case D 

8 

10 

0.05 

0.148 

1.00 X 10-2 

0.31 

5 


Note. — We give the number of protoplanets (iV), the surface density of protoplanets at 1 
AU (Si), and the semimajor axis of the innermost protoplanet (ai). The total mass of a system 
(Mtot), the dispersion of the eccentricity distribution ((e^)^/^), and the physical radius of the 
innermost planet {Ri/R(^) are given from above parameters. The last line of this table is the 
number of calculations in each iV-body case. 







Table 2. Initial conditions of three-planet simulations 



-^p/ 

Tfp/M0 

a-2 

[AU] 

b a2eini/ru 

da 

[AU] 

case 3A 

1 

1 

0.1 

4 

0 

5.04 X 10-3 

case 3B 

10V2 

1 

0.1 

4 

0 

5.04 X 10-3 

case 3C 

10 

1 

0.1 

4 

0 

5.04 X 10-3 

case 3D 

10-V3 

0.1 

0.1 

4 

0 

2.34 X 10-3 

case 3E 

10V3 

10 

0.1 

4 

0 

1.09 X 10-2 

case 3F 

1 

1 

0.13/2 

4 

0 

5.04 X 10-^ 

case 3G 

1 

1 

0 . 1 V 2 

4 

0 

1.59 X 10-3 

case 3H 

1 

1 

0.1 

5 

0 

6.30 X 10-3 

case 31 

1 

1 

0.1 

6 

0 

7.56 X 10-3 

case 3J 

1 

1 

0.1 

4 

2.5 

5.04 X 10-3 

case 3K 

1 

1 

0.1 

4 

5 

5.04 X 10-3 


Note. — The radii of planets the mass of planets (Mp/M,^), the semimajor axis 

of the middle planet ( 02 ), the orbital separations normalized by the Hill radius (b), the initial 
eccentricities (a 2 eini/rH) are given as initial conditions. The orbital separations between planets 
(da) is derived from da = brn = 6(2Mp/3M*)^/^a2. 
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Table 3. Results of iV-body Cases 



(n) 


(o-ll) 

[AU] 

(^11)/^esc 

{Mi2}/M^ 

{“i 2 > 

[AU] 

{^12)/^esc 

(^rem) / ^esc 

case A 

5.85 ± 1.01 

4.41 ±0.85 

0.20 ± 0.05 

0.18 ±0.08 

3.72 ±0.66 

0.17 ±0.06 

0.24 ± 0.26 

0.40 ± 0.32 

case B 

2.20 ±0.40 

16.7 ±2.8 

1.05 ±0.23 

0.16 ±0.14 

7.10 ±2.4 

1.04 ±0.56 

0.41 ± 0.35 

0.33 ± 0.0060 

case C 

3.20 ±0.98 

17.2 ±3.6 

0.11 ±0.03 

0.16 ±0.32 

11.7±3.4 

0.11 ± 0.05 

0.32 ± 0.20 

0.62 ± 0.53 

case D 

3.80 ±0.75 

0.0504 ± 0.0092 

0.060 ± 0.002 

0.47 ± 0.37 

0.0441 ± 0.0076 

0.057 ± 0.003 

0.36 ± 0.10 

0.44 ± 0.23 


Note. — Average values of the numbers of formed planets ((n)), the masses, semimajor 
axes, eccentricities normalized by their escape eccentricities of the largest planets {{Mu), {an), 
(eii)/eesc), and thoses of the second largest planets ((M/ 2 ), {ai 2 ), (ei 2 )/eesc) and the eccentricities 
normalized by their escape eccentricities of the other planets ((erem)/eesc) are denoted. 



